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Abstract 

The partition function of the g-state Potts model with random fer- 
romagnetic couplings in the large-g limit is generally dominated by the 
contribution of a single diagram of the high temperature expansion. 
Computing this dominant diagram amounts to minimizing a partic- 
ular submodular function. We provide a combinatorial optimization 
algorithm, the optimal cooperation algorithm, which works in poly- 
nomial time for any lattice. Practical implementation and the speed 
of the method is also discussed. 



Introduction 

Spin models with discrete symmetry are often used to describe order-disorder 
transitions, such as in magnetic systems or in absorbed monolayers, etc, for 
a review see in0. In two-dimensional regular lattices and with homogeneous 
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couplings there are some models (Ising model, critical Potts model, eight- 
vertex model, etc.) with an asymptotically exact solution for the partition 
function in the thermodynamic limit 0. In the presence of quenched, that 
is time independent, disorder this type of exact solutions are scarce which 
has greatly hampered our understanding about collective phenomena in dis- 
ordered systems, for example about the properties of spin glasses. 

For disordered systems, such as for spin models with random couplings 
and/or fields one usually performs an exact computation on finite lattices 
and the results are then extrapolated. An important limitation of this fi- 
nite lattice method is the available sizes, which can be treated by different 
computational methods. For most of the discrete spin models a computa- 
tion of the ground state represents a hard optimization problem when there 
is frustration and it is even much more complex to compute the partition 
function. 

In this respect the Ising model, which is equivalent to the Potts model 
with q = 2 states, is an exceptional problem since the ground (minimum 
energy) states can be found in polynomial time using the solution of max cut 
problems in planar graphs that involves matching theory H [[| ||. The full 
partition function of the Ising problem can also be computed for lattices of 
bounded genus ||. Another example of application of optimization theory 
to lattice statistical mechanics is the Random Field Ising Model where the 
ground state can be found in polynomial time in any dimension p. A 
discussion about optimal cuts in statistical mechanics can be found in || 
and a general review about the application of combinatorial optimization 
methods in statistical physics can be found in[l(| . 

For the q = 3-state Potts model, when the number of states is three, 
the computation of the minimum energy of the system involves (NP-)hard 
combinatorial problems such as minimizing the 3-cut in a graph, and no 
efficient procedure is known for computing the partition function. 

If the number of states can be arbitrary large (and given as part of the 
input) the partition function has an exponential number of terms, and the 
problem may look hopeless. In this paper we show that this is not the case: 
when the number of states tends to infinity, the partition function can be de- 
termined up to an integer multiplicative factor, and some average properties 
of the system can be inferred. This evaluation goes trough the minimization 
of a submodular function. Note that in [|ll|] a simulated annealing method 
and an approximate combinatorial method have been used. The solution we 
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give here is both exact and much more efficient. 

Even though there are ready-made algorithms for general submodular 
function minimization |T2| |T3| |TJ|] (T5J, the resources they require are far 
beyond the possibilities of the present application. While the existing first 
implementations of general submodular function minimization may require 
too much time already to solve the optimization problem on a lattice con- 
taining a dozen of sites (vertices of the input graph corresponding to the 
variables of the general problem), we want and will solve problems involving 
more than two hundred thousand 'sites'. For this we use methods inspired 
by |]16| ]17| [IB] |19j that fit the particularity of the application and develop 
an algorithm which is tailored to graphic submodular functions. Besides 
the computational advantages of this approach, we are also aiming at con- 
ceptual simplicity and elegance. Readers who seek for more mathematical 
background and connections can consulte pCfl . 

In Section [1] we describe the problem of physics and show how it reduces 
to a function maximization. In Section |2| we formulate the mathematical 
problem, the preliminaries and the broader context of our algorithm; in Sec- 
tion |3| and f| we show how to construct an optimal solution and in Section |5] 
we state the algorithm. Finally, in Section |6| we analyze the practical imple- 
mentation of the method and display a few computational results. The paper 
is self-contained, the algorithms and their correctness will be fully proved. 



1 The problem in terms of physics. 

For a review of the huge amount of work devoted to the Potts model see 
2jfl . In this section we recall the Potts model and show why its partition 
function (up to an integer multiplicative factor) is determined, if the number 
of 'states' tends to infinity by the minimum value of a submodular function. 

A lattice of spins is given at each site of which is a variable. We denote the 
number of sites by n. Each variable a\, . . . , a n can take values in {0, 1 . . . , q — 
1} for a given integer q. Pairs of neighboring sites on the lattice are called 
bonds, the number of bonds by m. 

Each configuration a = (o~i, . . . , cr n ) has an energy 

E{a) = Y. K iA^ (1) 
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where Oi G {0, . . . , q — 1} the sum runs over all bonds ij, Kij G R+ is a 
given non- negative weight of bond ij, and 5 a b is the Kronecker symbol (1 
if a = b and otherwise). For mathematicians: a lattice is a graph, sites 
are vertices, bonds are edges. The aim is to compute or approximate the 
partition function 

Z{K i5 : ij is a bond) = ^ exp(E(a)) = X) II e K " S ^ , (2) 

a a ij 

where the summation runs over all assignments of values a = (eri, . . . , a n ), 
and the product over all bonds ij (the inverse temperature has been included 
in the couplings K^). 

We follow 0. Note that e^ 5 = 1 + (e K - 1)5 for 5 G {0, 1} ; introduce 



Uij = e Kii — 1 and expand the product of sums: 

IK ' = 11(1 + (e*« - l)5 CTiCTj ) = X II 

where the summation runs over all subsets of bonds F. Substituting this to 
(D we get: 

o F ijGF 

=EEn^A, J =E? c(F) n^ (3) 

F a ijeF F ijeF 

where the sum runs over all subsets F of bonds, and where c(F) is the 
number of connected components of F on the set of all sites, counting also 
the isolated sites among the components. (By convention if F is empty then 
TiijeF v ij — 1-) We got the last equality by counting the number of different 
a = (o~i, . . . , cr n ) for which the product is nonzero; since it is nonzero if and 
only if a(i) = a (J) for every ij G F, that is, if and only if a is constant 
on every connected component of F, all possible a can be enumerated by 
choosing an element of Z q for every connected component of F independently. 
Therefore we have exactly q ^ such a. 

Note that this sum has 2 m terms, and that in (H) q does not need to 
have an integer value. (This provides a way of defining Potts' model for 
non- integer values [p^| .) Clearly, > is equivalent to > 0, and we 
can introduce a new set of variables a^- with 

Vij = q aij , 



4 



and the partition function becomes 

Z = Y,q c{F)+ ^ paiJ ■ (4) 

F 

Finally one introduces the function 

f(F) = c(F) + J2 ay (5) 
ijeF 

so that 

Z = Y,q f{F) - (6) 

F 

As pointed out in ref fLl| , while q tends to infinity the sum in @ is 
asymptotically equal to Z — Nq?* where /* = maxFc{i,..., n }f(F) and N is 
the number of optimal sets. Note that if the weights Wij are arbitrary reals, 
the degeneracy N is likely to be 1. 

We will analyze the combinatorial optimization problem of finding this 
maximum. The problem of finding the degeneracy of the optimal set will not 
be addressed. In the analysis of the Potts model magnetic properties are re- 
lated to geometrical properties of the optimal sets A, and specifically to their 
fractal dimension (if any). We assume that, on average and for large enough 
lattice, this fractal dimension is the same for all sets, and therefore finding 
a single optimal set is sufficient for our purpose. We will present methods 
which apply to any set of weights {a%j} in an arbitrary lattice (graph) and 
provide a simply stated and efficient algorithm. 



2 The problem in terms of graphs. 

For basic graph theoretic notions, terminology and notation we refer to . 
Given a graph G = (V, E) on n vertices with weight function w : E{G) — > R, 
we want to solve the following problem: 

Maximize {f G ,w( A ) = c g(A) + £ w(e) : A C E{G)} (POTTS), 

where cg{A) is the number of connected components of the graph G(A) = 
(V, A). Note a slight change of notation with Eq. |^ where an edge is denoted 
by ij the two corresponding adjacent vertices. We can suppose < w(e) < 1 
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for all e e E(G), because if w(e) = 0, then deleting the edge, if w(e) > 1, then 
contracting it (identifying its two endpoints) we get an equivalent problem. 

When no confusion is possible we will simply use / or f G for fo, w , c for 
cq. This function that has to be maximized in order to solve (POTTS) has 



a crucial and well-known property (see for instance [23|, Exercise 6.2). For 
any two subsets A and B of edges of G : 

f(A) + f(B) < f(A U B) + f(A n B). (SUPER) 

Indeed, it is clear that £ eeA w(e) +E e eB w(e) = EeeAuB w(e) +E ee AnB w(e), 
so we only have to show that c(A) + c(B) < c(A U B) + c(A n B)\ that is c 
satisfies the (SUPER) property. For the sake of completeness we include a 
proof. Proceed by induction on|A\S| + |S\A|. If this number is then the 
statment is obvious. If not, then there exists say e G A \ B. By induction 
c(A \ e) + c(B) < c(A U B \ e) + c(A fl B). Deleting an edge increases the 
number of connected components by at most 1 and so we are done unless 
the equality holds and c(A \ e) = c(A) and c(A U B \ e) = c(A U B) + 1. But 
this is impossible since c(A \ e) = c(A) means that the two endpoints of e 
are in the same connected component of G(A \ e), and then also in the same 
component of G(A U B \ e). 

A function g defined on the set of subsets of a set S, with values in R is 
said to be sub-modular if g(A) + g(B) > g(A U B) + g(A fl B) for all subsets 
A and B of S. The (SUPER) property shows that the function — /, whose 
minimization solves (POTTS), is submodular. 

By now several polynomial algorithms [12| |14| are known to 



minimize a submodular function, therefore we can already claim that there 
exists a polynomial algorithm solving our problem (POTTS). However, these 
general algorithms are not efficient enough for our instances. Fortunately, 
the specific properties of our problem make possible the use of a consid- 
erably simpler and quicker algorithm. Several more specialized algorithms 



were already known (see [20] for an extended presentation) and in particular 
Cunningham |I7| worked out a specialized algorithm to the "optimal attack 
problem" :"The weight of an edge represents the effort required by an at- 
tacker to destroy the edge, and the attacker derives a benefit for each new 
component created by destroying edges . . .". (POTTS) is equivalent to the 
"optimal attack problem" but with a 'complementary' viewpoint : the weight 
of an edge represents the benefits of cooperation between two vertices (say, 
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researchers) , furthermore there is a unit support for each component (say, for 
each research project). So there is a loss of support when two components 
unite. Cooperate optimally! 

The algorithm described below has its roots in [17j and is similar to [[IS 



f pTOf] . It has the best asymptotic worst case complexity, the same as and 
is probably the simplest, both conceptually and in the use of computational 
resources of the implementation. Sophisticated ingredients were necessary 
for finding this solution, but we don't need to make explicit use of them. 
The interested reader may find more details on that subject in |2(| . 



We return now to the (POTTS) problem. An important observation is 
the following: if A is a set of edges, and X is the vertex set of a connected 
component of G(A), then adding to A edges induced by X, that is with both 
endpoints in X, increases the value of /. Thus in an optimal solution A* 
of (POTTS) each connected component of A* contains all the edges of G it 
induces. Let us call Potts partition any partition of the vertices such that 
each class induces a connected subgraph of G. 

All we have to determine is : a Potts partition maximizing the total weight 
of edges induced by the classes plus the number of classes. 

As a consequence, following |Tj| , during all the procedure it will be suffi- 
cient to consider subsets of vertices instead of subsets of edges. 

A set F of edges optimizing (POTTS) will be called optimal, its associated 
Potts partition Vg{F) = {Xi, . . . Xk} will be called optimal as well. For any 
subset F of edges let w(F) = J2eeF w ( e ) an d for any subset X of vertices let 
w(X) = J2e€E(x) w ( e ) where E(X) is the set of edges with both extremities 
inX. 



3 Extension of a solution 

In this section we will see how any set of edges maximizing (POTTS) for a 
graph G minus a vertex, may be completed so that it maximizes (POTTS) 
for the graph G itself. 

Let G be a graph with a weight function w on the edges, let x be a vertex 
of G and let G' be the subgraph of G obtained by deleting x. The edges of 
G' keep their weight. Let F' C E{G') be an optimal solution of (POTTS) 
for G'. 
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The following Lemma is an easy consequence of the (SUPER) property. 

Lemma 1 There exists an optimal solution of (POTTS) for G which con- 
tains all edges of F' . 

Proof: Let F be any optimal solution of (POTTS) for G. By the (SU- 
PER) property one has in G: 

fa{F U F') > f G (F) + fc(F') - f G (F n F 1 ). 

Both F' and F n F' are subsets of edges of G' and since F' is optimal in 
G', f G (F') = 1 + f G ,{F') >\ + fo(FC\ F') = f G (F n F'), whence : 

f G (FUF')>f G (F). 

Since F is optimal in G by hypothesis, we get that F U F' is optimal too, 
and the Lemma is proved. □ 

As a consequence of Lemma [l] we can obtain an optimal solution F* for 
G by adding edges to F'. Let X±, . . . ,Xk be the vertices of the connected 
components of G'(F'), then V G (F') = {X 1} . . . , Xk, {%}}■ Each connected 
component in G(F*) will be either an element of V G (F') or the union of at 
least two elements of V G (F'). The next Lemma shows how the value of a 
Potts partition is affected when a subset of connected components are put 
together. 

Lemma 2 Let V and V' he two Potts partitions of G such that V = (V \ 
W) U {UXi; Xt E W} for some W C V. Then 

f G (V) = f G {V) - - 1 - w(E(W))), 

where E(W) denotes the set of edges of G joining vertices belonging to two 
different sets in W. 

In particular if V is optimal then |W| — 1 — w(E(W)) > and if V is 
optimal then |W| - 1 - w(E(W)) < 0. 

Proof: Replacing in V the sets of W by their union decreases the cardi- 
nality of V by | W | — 1, on the other hand let A and A' be the sets of edges 
induced by the classes respectively of V and V', A = A' U E(W) and so the 
weight of the edges augment by w(E(W)). The equality is then proved and 
the rest follows directly. □ 

We are now ready to prove the following theorem: 
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Theorem 1 For any W C Vg{F') containing {x} and minimizing |W| — 
1 - w(E(W)), the set F* = F' U E(W) is optimal for (POTTS) in G. 

Proof: By Lemma [l] we know that there exists an optimal solution F* 
of (POTTS) for G which contains F' . Let W C V G (F') such that W = 
Uxigw-^i is an element of Vg{F*). If W doesn't contain {x} then, since F' is 
optimal in G', we get by Lemma || that | W| — 1 — w(E(W)) > (this value is 
the same in G and G'). On the other hand F* is optimal in G which implies 
that \W\ - 1 - w(E(W)) = 0, but then F* \ £(W) is also optimal for C 
Hence there exists an optimal solution F* of (POTTS) for G containing F' 
and such that any element of Vg{F*) not containing x is already in V G {F'). 

So any W C Vg(F') containing {x}, such that W = Ux i ew- 5 Q induces 
a connected subgraph of G, and minimizing under these conditions |W| — 
1 — w(i?(VV)), will provide an optimal solution F' U E(W). (Notice that this 
minimum is < since for W = {{x}} we get 0.) 

In order to prove the theorem it remains only to show that we may skip the 
connectivity constraint : we show now that any subset W C Vg(F'), {x} G 
W, minimizing |W| — 1 — w(E(W)) corresponds to a connected subgraph 
in G. Assume not, and let Wi,...,Wj (/ > 2) be the subsets of V G (F') 
corresponding to the vertices of the connected components induced by W = 
Uxi^wXi in G, Wi being the one containing {x}. We have 

|W| - 1 - w(E(W)) = |Wi| + . . . + \Wi\ - 1 - io(£(Wi)) - ••• - w(E(Wi)) 

For i > 2, {x} Wi and so as already noticed |Wj| — w(E{Wi)) > 1, 
which implies \W\ - 1 - w(£(W)) > |Wi| — 1 — w{E{W x )) and so we get a 
contradiction to our assumption on W. □ 

At this point we see that any way to find W will provide a constructive 
algorithm for getting an optimal solution of (POTTS) : we start with a 
solution for a small subgraph (for example a one vertex subgraph) and then 
add the vertices one by one, computing at each step a new extended solution. 
Fortunately an optimal W can be computed and in the next section we are 
going to show how. 

We remark that the value |W| — 1 — w(E{W)) doesn't depend on the 
subgraphs induced by the subsets A; in W, so we may ignore these and 
work in a possibly smaller graph. To be precise : the result of shrinking 
the pairwise disjoint sets A 1; . . . , X k of Vg(F') in G is the graph shr(G) = 
(shr(V"), shr(_E)), where shr(^/) = {x, x%, . . . , x k } where X\, . . . , Xk are distinct 
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new vertices, and the function shr : V — > {x,x\, . . . ,Xk} is defined with 
shr(t> ) := Xi if v G Xi and shr(x) := x; the image of an edge e with extremities 
x and y is shr(e) with extremities shr(x)shr(y), moreover, if this is a loop 
(shr(x) = shr(y)), it is deleted; edges keep their weight, that is w s h r (shr(e)) = 
w(e) ; sets of vertices or of edges are replaced by the image sets. There is a one 
to one correspondence between the subsets W of Vg(F') containing {x} and 
subsets W of vertices of shr(G) containing x and for any such W = shr(W) 
one has | W| — 1 — w(E(W)) = \W\ — 1 — Wabx{W) (as defined in the preceding 
chapter w S h T (W) is the sum of weights of the edges in the subgraph of shr(G) 
induced by W). So W will be optimal if and only W = shr(W) minimizes 
\W\ - l-w shl (W). 



4 A network flow model 

Given a graph H with weight function w on the edges and W C V(H), let 
b(W) = \ W\ — 1 — w(W). From the previous chapter it is clear that any way 
of solving the following problem will provide an algorithm solving (POTTS) 
: given a vertex x in H, find W* a subset of vertices of H containing x 
such that b(W*) = min(b(W); W C V(H),x G W). This problem is solved 
both in Picard-Queyranne [E4| and Padberg-Wolsey EH], using a network 



flow model that will be described below. (The work pi |, p5[ solve much 
more complex problems involving several network flow computations - the 
subroutine below is the adaptation of an auxiliary procedure.) 

We first give some definitions and wellknown facts. Given a directed graph 
and a subset X of its vertices, S(X) will denote the set of arcs leaving X, and 
is called a cut; if s G X, t ^ X it is an (s,t)-cut; S(X, Y) denotes the set of 
arcs oriented from X to Y. A function c with nonnegative value on the arcs 
is called a capacity function. If F is a set of arcs then c(F) := J2eeF c ( e )- A 
network is a directed graph with capacity function. By a well known theorem 



of Ford and Fulkerson a minimum (s, t)-cut in a network may be found 



by computing a maximum flow in this network. Very efficient maximum flow 



algorithms are known |27 |. 



We now describe in several steps a network (D,c) = N(H,w) that will 
be associated to H and w: 

1. V(D) := V(H) U {s,t} where s, t are distinct new vertices; 
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2. Define for all u G V(H): 



p(u) := 1/2 w(uv), and then 

uey(H),«t)6£;(H) 

- if p(u) > 1 add an arc (s,u) of capacity c(s,u) = p(u) — 1, 

- if p(tt) < 1 then add an arc (u, t) of capacity c(w, t) — 1 — 

3. to each edge wt> of E(H) we associate the arcs (w, t>) and (t>, u) of capacities 
c(m, f) = c(v,u) = \w(uv). 

Let C V(H), 5({s} U W) is an (s,t)-cut of N(H,w) and reciprocally 
any (s,t)-cut of N(H,w) corresponds to a subset of vertices of H. There is 
a close relationship between the capacities of the (s, t)-cuts of N(H) and the 
values of b on subsets of vertices of H : 

c(5({s} U W)) = \W\ - w(W) + K = b(W) + K + 1, (CUT) 

where K := c(S(s)) = T,( s ,x)eAc(s,x). 

Indeed, let us see how the capacity of the cuts changes if we start with 
the set {s} inducing an (s,t)-cut of capacity c(S(s)) = K and then 'add' to 
it the vertices of W one by one: 

The contribution of adding v to the side of s is 1 — p(v), for either it 
decreases from p(v) — 1 to (this happens if p(v) > 1 see the first case of '2.' 
in the construction of N(H,w)), or it increases from to 1 —p(v). Thus the 
contribution of these arcs is: 

J2l-p(v) = \W\-w(W)-l/2 ]T w(xy). 

v&W xy&E{H),x&W,y£W 

On the other hand the contribution of the arcs between W and V(H) \ W 
is clear: at the beginning it is zero, and at the end it is: 

c(W,V(H) \ W) = 1/2 J2 w(xy). 

xy£E(H),x£W,y(£W 

The change comparing to K is provided by the sum of the two contribu- 
tions, which is \W\ — w(W), and so equation (CUT) has been proved. 
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From (CUT) we get that any (s, t)-cut C = {s}UW containing the special 
vertex x and of minimum capacity will correspond to our goal : a subset W 
of vertices of H containing x and minimizing b. Such a cut is easy to obtain : 
shrink {s, x} or equivalently, add to N(H, w) an arc of infinite capacity from 
s to x, and compute in the new network a minimum (s, t)-cut : this cut will 
contain x and will be minimum among (s, t)-cuts containing x. 

Let us note that this method can be easily generalized to minimize any 
modular shift of —w(W). Adding \W\ — 1 is just a particular choice of a 
modular function. That is, one can minimize (J2wew m ( w )) ~ w(W), where 
m : V — > R. Indeed, for generalizing from \W\ to m(W) one only has to 
write m(x) instead of 1 in both cases of Step 2 in the construction of N(H,w). 

5 The Algorithm 

In this section we state the 'Optimal Cooperation Algorithm' whose validity 
is a consequence of the results of the two preceding chapters. 

At each iteration a subset U C V and a partition V of U will be at hand. In 
Step we give trivial initial values; in Step 1, we choose an arbitrary vertex 
u to be added to U and through the following steps we compute a subset W 
ofVU {u} providing a partition of U U {u} consisting of the sets contained 
in V but not in W plus a new class containing all the other vertices ( that is : 
the new class is the union of all sets in W). The set of edges corresponding 
to this new partition is the union of the set of edges corresponding to the 
previous one and the edges between the sets in W. 

OPTIMAL COOPERATION ALGORITHM 

INPUT: A graph G = (V,E), and a weight function w: E — ► (0, 1), where 
(0,1) := {t E R : < t < 1}. 

OUTPUT: A partition V* optimizing (POTTS) : the set of edges A* induced 
by the sets in V* maximizes {/g,u>(-<4) := cq(A) + J2 e eA w ( e ) '■ A ^= E}- 

0. U :=<fr,V := 0. 

Do n times consecutively Steps 1 to 5, and then define V* = V: 

1. Choose a vertex u e V \ U. 

2. Define H and % as the result of shrinking the classes of V in G(U U {u}) 
with weight function w restricted to the edges of G(U U {u}). 
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3. Construct N(H,wh), add an arc of "infinite" capacity from s to u and 
using any maximum flow algorithm compute, in the so obtained network, a 
minimum (s, i)-cut defined by the set of vertices C = {s, u, x%, . . . , x k } where 
each Xi is a vertex of H corresponding to a set Xj G V. 

4. Define R := {u} U X 1 U . . . U X k 

5. Redefine U and V: U := U U {«}, V := (V \ {X u . . . , X fc }) U {R}; 
END; 

Remarks: 

- "Infinite" means big enough to avoid a minimum cut containing this 
arc: for example the sum of the capacities of all arcs. 

- Note that C can be equal to {s,u}. 

- The graphs H can also be constructed iteratively with only one shrinking 
in each iteration, by adding u and then shrinking the set C \ {s}. 

- The choice for u is completely free; this freedom could be used for 
making the computations easy. 

Theorem 2 The output V* is an optimal partition for (G,w). 

Proof: At step 0, V :— is an optimal solution for the subgraph of G 
induced by U := 0. Assume now that V is an optimal solution for the 
subgraph of G induced by U, and show that after applying steps 1 to 5 the 
new partition is optimal for U U {u}. From the preceding chapter we know 
that W := C \ {t} = {u, x%, . . . , Xk} is a subset of vertices of H containing 
u minimizing \W\ — wh{W) — 1. But H and Wh are obtained by shrinking 
the classes of V in (G(U U {«}), w), hence W = {u, X\, . . . , X k } is a subset 
of {u} U V minimizing |W| — 1 — w(E(W)) and by Theorem [I] the partition 
(V \ {X u . . . , X k }) U {R} is then optimal □ 

Notice that the "Optimal Cooperation algorithm" consists merely of n— 1 
network flow computations and n — 1 shrinking. 

6 Implementation and evaluation 

In this section we discuss some aspects of the practical implementation of the 
optimal cooperation algorithm. Here we consider a two dimensional square 
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lattice where the weights, w, follow a bimodal distribution: 



P(w) = p5(w — wi) + (1 — p)8(w — W2) (7) 

where 8(x) — 1 if x = and S(x) = otherwise. We restrict ourself to the 
symmetric case p = 1/2. As explained in [[Uj in the infinite lattice limit 
there is a second-order phase transition in the model at W\ + W2 — 1. In the 
strongly disordered phase, Wi + w 2 C 1 the optimal set of edge is A = 0, 
while in the strongly ordered phase, W\ + W2 1 it is A = E. In between, 
the optimal set A has a non-trivial structure and at the transition point it is 
a fractal reflecting the critical properties of the system [JTTJ . 



As explained in the previous section, for a graph with n vertices the 
algorithm consists in n — 1 steps. At each step a new site is incorporated 
to the current graph for which the optimization has already been performed. 
There is a freedom of the order of the added sites, which can be used to speed 
up the computation, as well as to obtain the optimal sets of intermediary 
lattices. For example from the computation of the optimal set for an L x L 
square lattice with free boundary conditions one can also obtain the optimal 
sets of L — 1 minors of size ranging from 1 to L — 1 with free boundary 
conditions, which is useful for finite size analysis. 

The "engine" of the algorithm is the max flow solver, which is called at 



each step. We have used the Goldberg and Tarjan algorithm pq| . Another 
important point concerns the building of the network (point 3 in the Optimal 
Cooperation Algorithm). It is possible to build it explicitly as another data 
structure than the lattice. Alternatively one can increase up to a sufficiently 
large value all the capacities of the edges of a spanning tree of each connected 
component as mentioned in the previous paragraph. In general the cpu time 
required to explicitly shrink the connected components is larger than the cpu 
time needed to simply increase some capacities. On the other hand, keeping 
all vertices will increase the size of the network in which the max flow has to 
be found. If all the capacities are large {w\ + W2 > 1, that is we are in the 
ordered phase) the network will have only two sites at each step, and therefore 
the effort to build the network will be largely compensated by the speed at 
which the maximum flow is found. In the opposite extreme, when the optimal 
set tends to be empty, it is questionable which of the two strategies is the 
best. We have found that explicitly shrinking the connected components is 
always the most efficient. As a matter of fact the computation time is about 
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L 


time 


16 


~ one second 


32 


~ half a minute 


64 


~ five minutes 


128 


~one hour and a half 


256 


~ one day 


512 


~ two weeks 



Table 1: Typical cpu time (see text) 



ten times larger at the critical point w\ + w 2 = 1 than at the ordered phase 
(wi + W2 > 1) and about three times larger than at the disordered phase 

(W! +W 2 < 1). 

As an illustration we display in the table below some typical cpu times 
needed to find an optimal set on a periodic square lattice of various size, L, 
at the critical point w± + w 2 = 1. From these results one can extract the 
empirical rule that in our algorithm the computation time grows approxi- 
mately as the fourth power of the linear size L, thus as the second power of 
the number of sites. These cpu times refer to a PC with a pentium III, 800 
megahertz processor. 

In conclusion in this paper we have presented and implemented a combi- 
natorial optimization algorithm which solves the optimal cooperation prob- 
lem in polynomial time. The procedure, which can be used to obtain the 
partition function of the g-state random bond Potts model in the large q- 
limit is, up to our knowledge, the first application of submodular function 
optimization in (statistical) physics. Detailed results about the g-state ran- 
dom bond Potts model, such as properties of phase transitions in two- and 
three-dimensional lattices will be presented in a separate publication [^|. 
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